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Abstract 

We derive a differential equation that is regular at the collision of two equal-mass bodies with 
attractive interaction in the relativistic action-at-a-distance electrodynamics. We use the energy 
constant related to the Poincare invariance of the theory to define finite variables with finite deriva- 
tives at the collision. The collision orbits are calculated numerically using the regular equation 
adapted in a self-consistent minimization method (a stable numerical method that chooses only 
nonrunaway solutions). This dynamical system appeared 100 years ago as an example of covariant 
time-symmetric two-body dynamics and aquired the status of electrodynamics in the 1940's by the 
works of Dirac, Wheeler and Feynman. We outline the method with an emphasis on the physics 
of this complex conservative dynamical system. 
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I. INTRODUCTION 



Delay equations play an essential part in Maxwell's electrodynamics because 
of the finite speed of propagation of the electromagnetic fields. Another more 
modern version of electrodynamics, action-at-a-distance electrodynamics [1], de- 
veloped in the 1940's as an alternative to avoid the divergencies of perturbative 
quantum electrodynamics [1, 2], was stopped by a main difficulty: delay equa- 
tions. The fundamental problem of a non-perturbative calculation of the level 
shifts of hydrogen depends on our ability to deal with delay equations. Given a 
Hamiltonian description for this electromagnetic two-body problem [3], knowl- 
edge of the orbits can be used in EBK quantization. Besides delay, an extra 
difficulty with the numerical calculation of an orbit for the relativistic two-body 
problem is the collision, where the equations of motion in usual form become 
singular. This obstacle has so far prevented the numerical study of this prob- 
lem. In this work we derive a delay differential equation that is regular at the 
collision of two equal-mass bodies with attractive interaction in the action-at- 
a-distance electrodynamics. The collision orbit is calculated numerically with 
a self-consistent minimization method that integrates the regularized equation 
[3] . We calculate numerically four collision orbits with energies from the nuclear 
to the atomic range. 

The usual Hamiltonian description of two-body dynamics is surprisingly restrictive within 
relativity physics: If Lorentz transformations are to be represented by canonical transfor- 
mations, only non-interacting two-body motion can be described. This is the content of 
the no-interaction theorem of 1964, which later in 1984 was proved for the local Lagrangian 
description as well [4]. A covariant version of Hamiltonian dynamics, constraint dynamics, 
was invented to overcome this group-theoretical obstacle, but it has a limited applicability 
[5], and in particular a constraint description of electrodynamics is not known at present. 
In the light of [4], the only available description of Lorentz-invariant two-body dynamics 
is via a Lagrangian built from the scalar invariant of the Poincare group, (the modulus of 
the separation four- vector) , that involves the time coordinates and yields delay equations of 
motion. This is a second, more modern confirmation that delay equations are in relativity 
physics and electrodynamics to stay. The idea to remove the field degrees of freedom from 
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electrodynamics goes back to Dirac [6] and later Wheeler and Feynman understood in 1945 
that act ion- at- a- distance electrodynamics was a theory that did not need renormalization. 
The subsequent program to quantize the two-body problem of the action-at-a-distance elec- 
trodynamics faced mathematical difficulties that can be summarized in one word: delay 
equations. History says that the famous seminar that never came from Wheeler (see Ref. 
[2], page 97) was due to difficulties in dealing with delay equations. In this same chapter 5, 
page 97 of reference [2], Feynman says that ' I didn 't solve it either — a quantum theory 
of half-advanced half-retarded potentials — and I worked on it for years... '. Since then this 
has been an open problem of atomic physics, with the main subsequent inputs coming from 
the works of applied mathematicians [7]. 

In 1903, Schwarzchild proposed a relativistic type of interaction between charges that 
was time reversible precisely because it involved retarded and advanced interactions sym- 
metrically [8]. The same model reappeared in the 1920s in the work of Tetrode and Fokker 
[9] and it finally became an interesting physical theory after Wheeler and Feynman showed 
that this direct-interaction theory can describe all the classical electromagnetic phenomena 
(i.e. the classical laws of Coulomb, Faraday, Ampere, and Biot-Savart) [1, 10]. Wheeler 
and Feynman also showed in 1945 that in the limit where the electron interacts with a com- 
pletely absorbing universe, the response of this universe to the electron's field is equivalent to 
the local Lorentz-Dirac self-interaction theory [6] without the need of mass renormalization 
[1, 11]. The Wheeler and Feynman program [2] to quantize the action-at-distance electro- 
dynamics and overcome the infinities of QED is still not implemented because of the lack of 
a Hamiltonian description [3]. As very little is known of this important physical problem at 
an analytical level, the knowledge of the trajectories can be useful in EBK quantization[3]. 

The isolated two-body system, away from the other charges of the universe is a conser- 
vative time-reversible dynamical system in the action-at-a-distance electrodynamics. We 
consider here only the equal- mass two-body system (mi = m 2 = m), henceforth called 1D- 
WF2B. The only postulate of the relativistic action-at-a-distance electrodynamics is that 
equations of motion be derived formally [12] by extremizing the parametrization-independent 
action 



S F = 



J mdsi — J mds 2 — e\e 2 J J 5(\ \xi — x 2 \ | 2 )ii • x 2 dsids 2 , (1) 

where Xj(sj) represents the four-position of particle i — 1,2 parametrized by its arc-length 
Si , double bars stand for the four-vector modulus — x 2 \\ 2 = (xi — x 2 ) • (xi — x 2 ), and 



the dot indicates the Minkowski scalar product of four-vectors with the metric tensor g^ v 
{doo = l,<7n = 922 = #33 = —1) (the speed of light is c = 1). The attractive problem is 
defined by Eq.(l) with e\ = — e 2 = e (positronium atom), while the repulsive two-electron 
problem is defined by Eq. (1) with e\ = e 2 = e. For the repulsive two-electron problem 
along symmetric orbits [— x 2 (t) = Xi(t) = x(t)], minimization of action (1) prescribes the 
following equation of motion 



where v(t) = dx/dt is the velocity of the first electron, of mass m and charge e, and r and 
q are the time-dependent delay and advance, respectively. The functions r(t) and s(t) are 
implicitly defined by the light-cone conditions 



In general, a neutral-delay equation such as Eqs. (2) and (3) requires a pair of world-line 
segments of trajectory as the initial condition (one world-line segment for each particle). As 
discussed in Ref. [13], the initial world-line segments can be provided in such a way that 
Eqs. (2) and (3) are well-posed, by using "maximal independent segments ". A pair of world- 
line segments is called independent if the end points of each segment lie on the forward 
and backward light-cones of a single point interior to the other segment. Last, a surprising 
existence theorem was proved for the symmetric motion of two electrons along a straight 
line [—X2(t) = Xi(t) = x(t)} (Eqs. (2) and (3) ). For this simple motion and for sufficiently 
low energies, it was shown in Ref. [7], that Newtonian initial conditions [ x(0) = x Q and 
v(0) = v ] determine the unique solution that is globally defined (i.e., that does not runaway 
at some point) [7]. Existence/uniqueness proofs are still lacking for the case of attractive 
interaction, and we hope that with the present regularization of the equations of motion 
such proofs can be facilitated. 

For the relativistic two-body system, the only known analytical solution is the circular 
orbit for the attractive problem [14, 15]. The first numerical method to solve Eqs. (2) and (3) 
was developed in [16] and converged to solutions up to v/c — 0.94. Later another method 
[17] converged up to v/c = 0.99. In reference [3] we developed a numerical method for 
the repulsive problem. Precisely because of the singularity, numerical methods and studies 





r(t) 



x(t)+x(t + q). 



x(t) + x(t — r), 



(3) 
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for the attractive problem are lacking, and the method developed in [18] has no hope of 
dealing with a near-collision. We hope that this work can start to fill this gap. This paper 
is organized as follows: In Section II we develop familiarity with the collision orbit and 
regularization issues. In Section III we study the behavior of the symmetric orbit near the 
collision with formal series expansions, and motivate the change of the evolution parameter 
(the time transformation). As this alone is not enough to accomplish manifest regularization, 
in Section IV we introduce the energy constants of the electromagnetic two-body problem 
to aid in the definition of two finite variables with manifestly finite derivatives (i.e. our 
regular differential equation). Because of the delay nature of the equations, coordinate 
transformations alone are not enough to prove that the derivatives are finite at the collision. 
In appendix A we make use of the energy constants to recognize the mathematical space for 
regular orbits and to prove the regularity of these derivatives. The material of appendix A 
provides an elegant alternative to the pedestrian construction of formal power series in the 
neighborhood of the collision (i.e., the material of section III). Last, in Section V we adapt 
the regular equation in a numerical method that integrates future and past histories until 
the eventual self-consistency of the histories. In this section we also calculate numerically 
several orbits in several energy ranges, and the implementation needed in each range to 
speed up convergence of the method. In section VI we put the conclusion and discussions. 

II. THE EQUATIONS OF MOTION 

Our regularization follows closely the Levi-Civita regularization of the Galilei-invariant 
Kepler problem [19-21], with the additional difficulties imposed by the Poincare invariance 
(i.e. delay). As with the Levi-Civita regularization, a time transformation alone does 
not accomplish regularization, and it is necessary to use the energy constant to remove 
infinities from the equation of motion. In the present Poincare-invariant case, besides a time 
transformation and use of the energy constant, it is further necessary to define special finite 
variables to accomplish manifest regularization. The non-local expression for the conserved 
energy of the electromagnetic two-body problem [1, 22, 23] is therefore used here in two ways: 

(i) As in the Levi-Civita regularization, to remove infinities from the equation of motion and 

(ii) to define the required finite variables with manifestly finite derivatives and to aid in the 
proof of regularity of these derivatives. Unlike the Levi-Civita regularization, because of the 
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delay nature of the equations, it is not possible to check the regularity just by performing 
coordinate transformations and taking limits, and the use of the energy expression provides 
an elegant way to perform these limits and to recognize the correct space of definition of the 
regular orbits, as discussed in appendix A. With that we accomplish manifest regularization. 

We henceforth use a unit system where m = c = e\ = — e 2 = 1. We assume that at 
t = —tc particle 1 is at x± — and moving to the right while particle 2 is at the same 
point and moving to the left (outgoing collision). The particles collide again at t — tc with 
ingoing velocities, so that x 1 (t) > 0, x 2a (t) and x 2b (t) < all along the unit cell of our orbit 
(see Figure 1). Because the transformations involved in the regularization are elaborate, 
we choose to work in the special Lorentz frame where the orbit is symmetric and therefore 
loose covariance in benefit of the intuitive picture. A covariant analysis shall be left for 
later work. In this work we consider only symmetric orbits of the equal-mass attractive 
1D-WF2B, whose equation of motion for particle 1 is 

dVi 1 f 1 {I -Via) 1 {l + V2b) \ (, 2X3/2 

dh 2 \rl(l+v 2a ) + rl(l-v 2b )j {i V±) ' [V 
while the condition of symmetric colinear motion defines the trajectory of particles 2 as 

x 2 (t) = - Xl (t) . (5) 

In Eq.(4) Vi stands for the instantaneous velocity of particle 1 (present time is the time t\ 
of particle 1), while v 2 b and v 2a stand for the velocities of particle 2 at the retarded and 
advanced light-cones, respectively. As with Eqs.(2) and (3), the electrodynamic interaction 
in Eq. (4) connects points that are in light-cone condition, as defined by 

r a = \xi (ti) - x 2 (t 2a )\ = c(t 2a ~ h) , (6) 
r b = \xi {h) - x 2 (t 2b ) | = c (f i - t 2b ) . 

In all the above, subscripts a and b indicate future and past times of particle 2 in light-cone 
to the present time t\ of particle 1 (see Figure 2). By use of Eq. (6), we can also express t 2a 
and t 2b as 

tia = h +r a /c, (7) 
t2b = h - r b /c, 
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and from Eqs.(6) and (7) we can derive the equation of motion for t 2a 

dt 2a (l + ^i) 



(8) 



dti (l + v 2a ) 

In the same way that we derived Eq. (8), the motion of x 2a , x 2 b, and t 2 b can also be derived 
from Eqs. (6) and (7) as 

dx 2a (l + ui) , n x 

dti (1 + v 2a ) 

dx 2b (1 - vi) 

~ v 2b - 



dh (1 - v 2b ) ' 

dt 2b (1 - ^i) 
dh {l-v 2b ) , 

Eqs. (4) and (8)- (9), together with the definitions of r a and r b of Eq. (6) constitute the 
complete delay equation that we consider in this paper. Along a symmetric orbit ( as defined 
by Eq. (5)), the retarded and advanced velocities v 2b and v 2a are defined from the velocity 
Vi(t) of particle 1 by 



v 2a (t) = - Vl (t + r a /c), (10) 
v 2b (t) ee -v 1 (t -r b /c) , 

as illustrated in Figure 2. By use of Eqs. (5) and (6), we can also show that along a symmetric 
orbit the two quantities r a (ti) and r b (ti) are defined by a single function r (t) as 

Ta if) =r(t), (11) 

r b (t + r a /c) = r(t), 

which is illustrated in Fig. 3. As the force is always attractive from both retarded and 
advanced positions in Eq. (4), after the velocities have switched opposite the interparticle 
distance must approach zero until the collision happens. One could conjecture of orbits 
where the particles reach the speed of light even before the collision. Such orbits, if they 
exist, will not be studied here. We shall show below that both velocities must tend to the 
speed of light as the particles approach the collision. 

To prove that the velocities must go to the speed of light at the collision, we need some 
monotonicity properties: From Eq. (4) it follows that ^ < for all times, and by symmetry 
one also has > and > , such that we can establish that 

al2a dt2b 1 

\v 2a \ <v 1 < \v 2b \ , (12) 
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for any time — tc < h < tc- It is also easy to show that r a (ti) and r b (ti) are piecewise 
monotonic functions of t\. From Eqs. (6) and (9) it follows that 



dr a = vi - v 2a , , 

dh ~ 1 + v 2a ' 1 J 



(14) 



and 

dr b _ vi - v 2b 

dti 1 — V 2 b 

As the velocities are globally monotonic, there is a maximum radius r (see Figure 1) attained 
at ti = t > —t c ,when Vi (ti = t ) = v 2a (ha), (such that ^ = 0). It must be that 

dv 

dh > ° 5 tl<tQ ' (15) 
dr a 

< ; h>t . 

dt 

As the collision happens at t\ = —tc < to, we can restrict to the increasing part of r a , 
(—tc < t < to) ) an interval where Eq. (11) determines the bound 

r b (t)=r a (t-r a /c) <r a (t). (16) 

For the complementary interval before the next collision (t c > t > t ), we have r a < r b . 
Because of inequality (16), when the largest radius r a (t) goes to zero, r b (t) must go to zero 
as well, such that this largest radius becomes the natural control parameter for the dynamics 
in the neighborhood of the collision. Our next Lemma shows that velocities v 2a ,vi, and v 2b 
must all tend to the speed of light in modulus when this largest radius r a goes to zero. 
Lemma: 

Assuming that a continuous solution (r (t) ,v (t)) exists in an open neighborhood of the 
collision point r a = and t = —tc , then we must have that both velocities go to the speed 
of light at t = —t c 

Proof by contradiction: 

Dividing Eq. (4) by Eq. (14) we obtain the following equation for the evolution of vi 
with r b 

dVl _ (l~V 2b ) / (l-V 2a ) | (l+V 2b ) \ ^ ^3/2 ^ 



dr b 2 (vi - v 2b ) \ r 2 a (1 + v 2a ) r\ (1 - v 2b ) 
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If neither of the velocities goes to the speed of light at the collision, then by Eq. (10) it 
must be that 

Vl = v c < 1, (18) 

Via = ~V C , 
V2b = -V C . 

Given that the velocities are bounded as of Eq. (18) and because of Eq.(16), the second 
term on the right-hand side of Eq. (17) dominates and integration of this dominant term 
yields 

v 1 = -+ (19) 

n 

with k a nonzero constant. Equation (19) predicts that V\ becomes infinite as r& goes to 
zero, an absurd, as we have assumed that v\ — v± < 1. We conclude that v\ (r a = 0) = 1 
and v 2 b (n = 0) = v 2a (r a = 0) = -1 . 

III. THE TIME TRANSFORMATION 



To motivate our regularizing time-transformation, we assume that at least one regular 
orbit exists and construct its formal series expansion in the neighborhood of the collision. 
For this we develop the function 1 + v 2a in a power series of r a and the function 1 + v 2 b in a 
power series of r b 

l + v 2a ~ ar q a + ai r q a +1 + ... 

l+v 2b ~ Br s b + B ir l +1 + ... (20) 

where s and q must be positive because the velocities have a bounded modulus (|t> 2 | < 1 ) 
and a and B > are to be determined later. It is easy to verify that the evolution parameter 
u defined as 

dh =du(l + v 2a ) (1 - v 2h ) , (21) 

regularizes Eqs. (8)-(9) and Eqs. (13)-(14) at the collision. The only problematic regu- 
larization left is the right-hand side of Eq. (4), which involves two indefinite limits at the 
collision. In the following we show that the second term on the right-hand side of Eq. (4) 
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vanishes at the collision. For this we divide Eqs. (4) and (14) by Eq. (13) and obtain 
differential equations for v\ and r b in terms of the evolution parameter r a 

**- (1 -" f)3/2 ^i-, 2a) + (1 + ?' (1 + ^ )r n (22) 



dr a 2 r\ {y x - v 2a ) \ '" r l{ l ~ ^2b) 

dr a K l-v 2b ,{ v 1 -vJ 2 + V 2 J 2 a + -- 1 6) 

Eqs. (22) and (23) are converted into ordinary differential equations by use of the series of 
Eqs. (20). Integrating the leading term of Eq. (23) we find 

r b ~ a J a - (24) 
2(g+l) 1 ; 

To obtain a series for V\ from Eq. (22), we start by noticing that the first term on the 
right-hand side of Eq. (22) is approximately 2 in the neighborhood of r a = 

(1 - v 2a ) ~ 2. (25) 

By use of Eqs. (20) and (24) we can show that the second term of the right-hand side of 
Eq. (22) is proportional to r^ s+s ~ q 

(l+^Hi+^j a rl s +s - q a (i + ^) (26) 

r\ (l-v 2b ) (l+v 2b ) 

As the velocity is monotonic, the limit of the ratio i n Eq. (26) must either be zero or at 

the worst this limit can be a constant value, which implies that qs+s—q > 0. In the following 
we assume that Eq.(26) vanishes near the collision, and obtain a leading approximation to 
the equations of motion. This leading approximation in turn calculates qs + s — q = 3, 
showing that the assumption is consistent. The pathological option qs + s — q = must be 
analyzed separately, as then Eq. (26) has a finite limit. This analysis again determines that 
qs + s — q = 3, showing that a vanishing limit for Eq. (26) is the only consistent choice. 

With the above in mind, the leading terms of Eqs. (4), (13) and (14) in the neighborhood 
of u — are 

t - 4 - < 27 > 

dvi 4\^2 ^ 3 /2 



du 



a 



dr b 

^~2(l + . 2 „). 
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Choosing u — at t — —tc and using Eq. (20) to eliminate v 2a , Eqs.(27) can be integrated, 
yielding 

r a ~ Au, 
1 + v 2a ~ au 9 
1-vi ~ 32m 2 , 

Eqs. (28) and (20) predict also the leading dependence of v 2b with u 

l+v 2b ~Brl = ^^u^\ (29) 

The last consistency condition on the solution is that x 2a (t) and x 2b (t) must describe the 
same orbit. This is accomplished if there exists a shift function Au (u) > 0, such that 

tl(u + Au(«)) = *2o(«), 

1 - Di (u + Au (u)) = l+v 2a (u) , 

l+v 2b (u + Au(u)) = 1- Vl (u). (30) 

Using the approximations of Eq. (28) and (29) we can solve Eq. (30) at the leading order 
of approximation with a term for Au (u) given by 

A M («) = 2 1/ y, (31) 

and also B — s — 2, q — 1 and a = 2y/2. Finally we can express the velocities in terms of 
the radii using Eqs. (28) and (29) 

l + v 2a ~ 2V2r a + ... , 
l-ui ~ 2r 2 a + ..., 
1 -ui ~ 2v^r 6 + ..., 

l+w 26 ~ 2r 2 + .... (32) 

For later use, it is interesting to obtain a further term of the series for v±, (Eq. (32)) in 
the following way: By use of the leading terms of Eqs. (28), we can express Eq. (22) as 

dr a ~ rUl + v 1 ) +U[ra) - {66) 
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The solution of Eq. (33) with the condition that v\ — 1 at r a — is 

-w-^IfeH^ (34) 

where C is an integration constant to be determined later. Expanding Eq. (34) in powers 
of r a we obtain 

l- Vl ~2r 2 a + 4Cr 3 a + ..., (35) 

which exhibits the next term in the series for v\ of Eq. (32). Eq. (35) also prescribes the 
next term in the series for v 2b of Eq. (32), which can be obtained by replacing r a by r b and 
V\ by — 1>2& in Eq. (35) 

1 + v 2b ~ 2r b 2 + 4CV b 3 . (36) 

This symmetry along time-reversible orbits is illustrated in Figure 3 (i.e., that v 2b (r b ) = 
—vi{r a — r b ) ). The constant C is related to the energy E of the orbit and is calculated 
in the next section. Last, it is of interest to notice that the term we disregarded in the 
approximation below Eq. (25), corresponding to the information from the past in Eq. (22), 
contributes to the expansion of V\ only at 5th order in r a . 



IV. TWO FINITE VARIABLES DEFINED BY THE ENERGY 

The conserved energy of the Kepler problem is simple and well known, while the corre- 
sponding energy for our relativistic problem is somewhat unfamiliar. The Poincare invari- 
ance of the Fokker Lagrangian determines a four-vector constant of motion, which involves 
an integral over a light-cone of the orbit [1, 22, 23]. For the one-dimensional symmetric mo- 
tion of equal masses, as explained in Ref . [23] , the total energy E^ F = E 1 + E 2 = 2E can 
be simplified in two independent constants (a time-reversed pair, Ei — E 2 — _£/) 

El = 7T^"^"(rr^)' (37) 
E2 = 7r^f"i + (T^)' (38) 



with Y a and Y b given by 
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where r+, v 2 and r_, v 2 stand for the radius and velocity at the advanced and retarded light- 
cones of particle 1 respectively, as illustrated in Figure 4. The total linear momentum of a 
symmetric orbit can also be expressed as P^f F = Pi + P 2 = [23] , with the time- reversed 
pair of constants P\ and P 2 given by 

p _P ^ 1 1 dt (l+ Vl _.) 1 rft (l+m-) 

2 ~ 1 2r a 2 7 rl (1 - Wl _) 2 /, r 2 (1 - Wl _) ' ^ 

where v 2 (ti) = —v\ (ti) and ff and t> f stand for the velocity at the retarded and advanced 
light-cones respectively. Using the equation of motion (Eq. (4)), one can check that the 
time derivative of either Eq.(37) or Eq. (41) vanishes. This implies that ^ = ^ = 
along the motion. The same applies to Eq. (38), which is the time-reversed of Eq.(37), and 
also to Eq. (42), the time-reversed of Eq. (41). Notice that the energy Eqs.(37) and ( 38) 
have the correct Coulombian limit far away from the collision, as Y a and Y b vanish for large 
distances. Using the near-collision behavior of Sec. Ill one finds that Y a and Y b as defined 
by Eqs. (39) and (40) have finite limits at the collision, which is the main motivation to 
introduce these Y variables as above. Recalling that Eq. (4) is the only equation left that is 
not regularized simply by the time transformation (21), in the following we use the energy 
Eqs. (37) and ( 38) to obtain regular differential equations to replace Eq. (4). 

The value of Y b as defined by Eq. (40) can be approximated by elimination of dt using 
the dominant term of Eq. (4) 

\r r 7^v- + °«)- < 43 > 

1 Jt 2b r + (J- + v 2 ) J- V2b (1 - v f) ' 

Evaluating the integral on the right-hand side of Eq. (43) yields an approximation for Y b 
near the collision 

n^-^ + ^=4, (44) 

a finite expression by use of Eq. (32). Substituting this expression into the energy Eq. (38) 
yields the approximation 

E ~ ~ + _L_ , (45) 



2r b y/l-v 



which predicts the same behavior for v 2b (r b ) as Eq.(32). Also the near collision behavior of 
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Y a can be derived from Eq. (32) 

Y (l+vi)(r a -r b ) (1 + Vl ) y/l- Vl 

2 2^2 K } 

At this point it is interesting to reverse the above argument; noticing that we could have 

derived the series for v\ ( Eq.(32)) from Eq. (44) with the assumption that Y b is finite. In 

Sec. Ill we had to use the existence of the shift function (Eq. (31)) to accomplish this same 

result. The fact that Eq. (32) can be obtained from the more elegant hypothesis that Y b is 

finite all along the orbit suggests the natural mathematical space for a regular solution; the 

set of orbits with finite values for Y a and Y b . Using Eq. (38) and the approximations of 

Eqs. (32), the constant C entering into Eq. (34) can be calculated as 

C = -2E. (47) 

The evolution of vi with respect to the parameter u can be calculated from Eq. (4) 

where £ is given by 

Z= (1 - vl) ( + + (4Q) 

V r a T 'b J 

It can be shown by use of Eqs. (32) that £ as defined by Eq. (49) has a finite limit at 
the collision. Notice the two numerically prohibitive features with Eqs. (48) and (49) : (i) 
Eq. (49) has indeterminacies of type , and (ii) Eq. (48) does not evolve starting from the 
initial condition v± — 1. An alternative way to obtain v\ is to solve the energy Eqs. (37) and 
(38) for v\. As the energy equations involve the variables Y a and Y b , it must be integrated 
along with equations for Y a and Y b . The needed equations of motion for Y a and Y b following 
from Eqs. (39) and (40) are 



dY a (vi - v 2a ) (1 - v 2b ) (1 + vi) _i( vi _ fl - vA 
du - 2rl 2 V(l-^i) a \ l+vj> [M) 

dY b (vi - v 2b ) (1 + v 2a ) (1 - v-i) £f v x 11 + vA 

*7 = 2V1 + H(TT^) +n VT^J- (51) 

The system formed by Eqs. (50) and (51) still contains indeterminacies of type jj, so we are 
not done yet. In Appendix A we show that the right-hand sides of Eqs. (50) and (51) have 
finite limits at the collision if: (i) the orbit has a finite energy and (ii) the variables Y a and 
Y b as defined by Eqs. (39) and (40) are finite along the collision at t — —t c . 
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V. NUMERICAL CALCULATIONS WITH THE REGULARIZED EQUATIONS 



We use here a numerical method previously developed by us in Ref.[3] for the repulsive 
case. This method approximates v 2a and v 2 b with two power series and integrates the 
regular equations (60) and (63) together with the regularized versions of Eqs. (8)-(9). This 
predicts future and past histories for particle 2, x 2a {t 2 ) an d x 2 b(t 2 ). Next a minimization 
scheme modifies the approximation for v 2a and v 2 b to improve the consistency of the two 
histories, until the eventual convergence to a consistent history for particle 2 is reached. 
Some observations are in order : (i) In Appendix A we show that the right-hand sides of 
Eqs. (50) and (51) have finite limits along a symmetric non-runaway orbit, henceforth called 
T/ and T/, respectively. The approximation for v 2a near the collision must be postulated 
such that that these limits are satisfied (which is expressed by Eq. (68) of Appendix A ). (ii) 
To describe a symmetric orbit in the most economical way, the global time-reversal symmetry 
can be embedded in the approximation for the velocities. According to this time-reversal 
symmetry, the advanced and retarded velocities must satisfy 



To satisfy (i) and (ii) and the near-collision behavior of Eqs. (32) we approximate v 2b and 
v 2a with an arbitrary function 9(r a ,r b ) such that 



If 0(r a ,r b ) is regular and evaluates to one at the collision, the above ansatz of Eq. (53) 
guarantees that T 6 of Eq.(64) is explicitly finite (i.e. there is no division of zero by zero to 
spoil the numerical calculations). Notice that under time reversal r a and r b are interchanged 

ij'a < ► ?"&) and V\ is exchanged by —v\ (i>i — > — V\), such that the two lines of Eq. 

(53) are exchanged. This is the embedding of the time-reversal symmetry v 2a < — > —v 2 b- 
For numerical convenience, the function 0(r a ,r b ) must be postulated in two different ways 
depending on the energy of the orbit: 



v 2a (t) 



v 2h (-*) , 



(52) 



r a (t) = n (-t) . 




(53) 
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(I) For orbits of E « 1, it is convenient to use a rational Pade approximation defined 

by 



0(r a , n) . -riT^^r (54) 



l + K l r b + K 2 rl + ... + K N r £ 
1 + k x r a + k 2 r 2 a + ...k N r% 
Eq.(54) is a quotient of a polynomial on r b over a polynomial on r a , which is constructed for 

the following reasons: (a) On the first collision, at t = —tc, the Tl (as defined by Eq.(64) in 

Appendix A) is explicitly finite, (b) On the second collision, at t — tc, the regularization as 

built with Eq. (54) is automatic because the embedded time-reversal symmetry exchanges 

r a and r b , such that T/ (as defined by Eq. (62) in Appendix A), is explicitly finite in the 

same way. 

(II) The atomic energy range, of greatest interest to physics, has E £3 1 , which is a 
difficult limit for approximation (I). Most of the counter-intuitive features of shallow energy 
orbits happen because the function 6 jumps abruptly from the value of one to another 
constant value at the Coulombian limit (r a ~ r b ^> 1). The numerically correct procedure for 
shallow energies is to postulate 6{v) with Spline interpolation [25] on the interval — 1 < v\ < 1 
, defined such that 9 — 1 at the collision (v — 1) and such that T[ and T b are finite at both 
collisions. This is accomplished with the definition 

6 = l + Vl-v 2 P(v). (55) 

We use up to 22 intervals to approximate the function P(v) with Splines[25]. Either the 
coefficients hi, /cat and K±, of the Pade approximation or the polynomial coefficients 
q of the cubic Splines are to be determined by the self-consistent minimization in each case. 
After these coefficients are substituted into Eqs. (53), and then into Eqs. (60) and (63) our 
regular equations of motion become ordinary differential equations that are integrated with 
a standard 9/8 explicit Runge-Kutta pair, generating the future and past of particle 2. Our 
self-consistent method calculates two functions that should vanish along a symmetric orbit 



Si(M) = x 1 (t)+x 2a (t), (56) 
S 2 (k, t) = x x {t) + x 2b (t) , 

at about m points along the orbit (m ~ 400). The interpolation coefficients (either k i: Ki or 
Cj) are changed by a least-square minimization algorithm (Levenberg-Marquardt algorithm) 
[26]. Notice that if we could find an analytical solution for the orbit, these coefficients could 
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be calculated by setting Sip to zero. In practice we determine a numerical zero for Sip of 
size 1 x 10~ 5 (see Table 1 ). As discussed below Eq. (36), only at O(r^) the information of 
the past becomes important near the collision, such that we should use N > 5, to include 
past information into Eq. (22). 

We are now ready to start the integrations from r a = r b = 0, using the approximate 
coefficients and E already calculated. The complete set of equations includes Eqs. (51) 
and (50) together with Eqs. (8) and (9) in terms of r a , r b) Y a and Y b . The initial condition 
r a = 0,r b = 0,v = 1, expressed in terms of the regular variables reads 

Y a (u = 0) = 0, (57) 
Y b (u = 0) = y/2. 

The energy E is a parameter that appears explicitly in the regular equation of motion, 
and the numerical procedure fixes E while the interpolation coefficients are adjusted by the 
minimization scheme. The velocity v± is calculated numerically by solving Eqs. (37) and 
(38) for v-\_. The numerically calculated orbits using the Pade approximation of Eq. (54) are 
shown in Figure 5 for four different energies ( E = —1.0 , E = 0.1, E = 0.5 and E = 0.8 ). In 
Table 1 we list the quantities related to these numerically calculated orbits. Notice in Table 
1 that the energy becomes negative when the maximum light-cone radius is lesser than the 
classical electronic radius. The orbits at atomic energies (E < 1) have two clearly separated 
regions: (i) a near-collision region where the velocity is very close to one and (ii) the turning 
region, where r a ~ r b ^> 1 and 6 — > ^- . This last segment of such orbits approximates the 
turning region of Coulombian orbits. In the collision region a relativistic orbit must deviate 
from a Coulombian orbit, because v i — > oo on the collision for Coulombian orbits. The 
transition between these two regions is abrupt, as illustrated in Figure 6. In Figure 7 we 
magnify the region of discontinuity for various energies, illustrating that the discontinuity 
in v changes shape with increasing E. In Figure 8 we show the numerically calculated orbits 
using the Spline interpolation for 9(v) ( Eq. (55)). The energy E and the relative error of 
these atomic orbits are shown in Table 2. Notice in Table 2 that the value of rd is converging 
to 3/2 as it should for Coulombian orbits. 



VI. CONCLUSIONS AND DISCUSSION 
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We derived a differential equation that is regular along the collision of two equal masses 
with attractive interaction of the action-at-a-distance electrodynamics, allowing the nu- 
merical study of these orbits for the first time. Our regular numerical method starts the 
integration exactly from the collision. Our procedure is not covariant because we restricted 
the work to symmetric orbits, having the energy as the only free parameter (the other pa- 
rameter should be the Lorentz boost parameter). A covariant treatment shall be left for 
future work, along the lines of Appendix B of Ref . [3] . The numerical results of Ref . [3] sug- 
gest that at least for the repulsive equal-mass case a boosted symmetric orbit is already the 
general nonrunaway solution. The generalization of this numerical study to the attractive 
case awaits a covariant regularization. 

Some failed attempts taught us that the different-mass attractive case is much more 
involved, and possibly not even regularizable. At present we do not even know how to coin 
a formal series solution near the collision for this case. One reason for that is the complexity 
of the energy expressions analogous to Eqs.(37) and ( 38) [23]. This fact confuses the 
definition of the finite variables for the different-mass case and the problem needs further 
study. 

The two-body problem with repulsive interaction also displays a singularity at high en- 
ergies, but surprisingly enough this is not because the particles collide (they never do). The 
singularity appears in the right-hand side of Eq. (2) when the particles come to the speed of 
light ( the denominator containing (1 — v) vanishes). It was found by several authors, first 
in Ref. [27] and later in Ref. [17], that the particles reach the speed of light and turn back 
keeping a minimum distance of approximation of about one classical electronic radius. The 
regularization of this problem shall be published elsewhere [28]. This minimum distance 
of closest approach of about one classical electronic radius is a kind of exclusion principle 
of the action-at-a-distance theory. It is of interest to notice that this exclusion behavior is 
already found with the post-Galilean low-velocity approximation to the Fokker Lagrangian, 
the Darwin Lagrangian [22]; The algebraic-differential equations of motion of the Darwin 
Lagrangian were studied analytically in [29] and the phenomenon of closest approximation 
was discovered. The distance of closest approximation with the Darwin Lagrangian is found 
to be exactly the classical electronic radius, as well as with the Fokker Lagrangian [27]. 
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VII. APPENDIX A : REGULARITY OF THE Y-DERIVATIVES 



In the following we show that the right-hand sides of Eqs. (50) and (51) both have a 
finite limit at the collision if: (i) the orbit has a finite energy and (ii) variables Y a and Y b as 
defined by Eqs. (39) and (40) are finite along the collision at t = —t c . It is convenient to 
introduce the variables z = y/1 — v\ and w = y/1 + v i and to define two manifestly finite 
quantities a and (5 as 

1 +v 2a 



a 



\/l - Vi 



Notice that z = and w = \/2 at v\ — 1, and the finite limiting values of a and (3 at the 
collision as predicted by Eqs. (32) are : a = 2 and f3 = 1/4. A concise way to obtain the 
expansions of (32) is by solving the energy Eq. (37) for r a and the energy Eq. (38 ) for r b 

Ta 2 (-w 2 zE -zY a + w)~ 2 Z) 

= 2 

Th 2 {-z 2 wE + wY b + z) ~ 2Y b Z ' 1 ' 

The variables p a and p b as defined by Eq. (59) have finite limiting values at the collision, 

p Q a — 1 and p° = 1. It is nice to observe that these finite limits are a consequence of E , 

Y a and Y b being finite in (59). We arrive then at a concise definition of a regular orbit: one 

defined by finite values of E , Y a and Y b . The derivative of Y a (Eq. (50) ) can be expressed 

in terms of these finite quantities as 



%±=Tl + z U- 2 + ^Y a + 2w*a(3^) - + (60) 

du \ pi 2w p{) p\ 

+z 3 ( 2c*4 - 2w 2 a^) + 6^ 4 4 - ^ 5 a^r, 

V pI pi) pI p\ 

where £, defined in Eq. (49), has the explicitly finite expression 

t _ ( I6(l-az/2)(l-f3z*) AY b a(3z^w 2 \ 

4 V pI pi J' [ } 

At z = the only nonzero term on the right-hand side of Eq. (60) is 



T[ = 4(1- z(wE-Y a /w))\ 
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(62) 



which is finite at the collision (z — 0, w — 0), such that the derivative of Y a as written in Eq. 
(60) is manifestly finite at the collision. The Y b derivative can be obtained by substitution 
of Eq. (59) into Eq. (51), and some algebraic manipulations yield 

^ = T( + (16 + + 8wEaY b 2 - Awa^ + -aY a Y 2 ) + (63) 
aw \ 2w Pa w J 

+z (-^- - 8EaY a Y h 2 - lQwE + 2a\ - Aw 2 E 2 aY 2 - A,aY 2 Y 2 
\wY b pi w 2 

( Y 2 El 

+z 2 2w 3 a[3^ T - 16- - — e 



Yb pI pI J pI 



where 



T( = ^(2w- aY b ) . (64) 

z 

The term Tjf defined by Eq. (64) was singled out in Eq. (63) because it contains a division 
of zero by zero at the collision, and is the only part of Eq. (63) that could in principle be 
singular. To show that T/ has a finite limit along an orbit of finite energy, we first use the 
energy Eq. (38) to express z as 

(r b + yWgg + r b + 1w*EYm) 
Z ~ w(l + 2Er b ) • (65) 

Substituting the approximation of Eq. (44) for Y h 

Y b ~V2-- + 0(z 4 ), (66) 

into Eq. (65), we find the expansion for z 2 (r b ) = (1 — vi) to be 

z 2 = 2r b V2 - AV2Er 2 b + O (r 3 b ) (67) 

Last, as discussed above Eq. (36) and illustrated in Figure 3, the particle exchange symmetry 
along symmetric orbits implies that the function Vi(r b ) is equal to the function — v 2b (r a ). 
Changing the argument of Eq. (67) from r b to r a and eliminating this r a in favor of V\ with 
Eq. (32) we obtain 

a = i±^£ = 2 - 2V2E z + O (z 2 ) (68) 

Substituting Eq. (68) for a into the formula for T/ , Eq. (64), and expanding we can 
determine that T b is finite and nonzero at the collision. 
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IX. FIGURE CAPTIONS 

Fig. 1 A symmetric trajectory, arbitrary units. Indicated are the maximum light-cone 
distance r , the velocity v\ of particle 1 at time t\ and the corresponding retarded and 
advanced velocities of particle 2, v 2b and v 2a . The trajectories of particle 1 (solid line on 
the right-hand side) and particle 2 (solid line on the left-hand side) are illustrated from the 
outgoing collision at t — —tc until the ingoing collision at t — tc- 

Fig. 2 The orbits of particle 1 (solid line on the right-hand side) and particle 2 (solid line 
on the left-hand side) in arbitrary units. Indicated is the symmetric point v\ — at t — 0. 
Also shown are the the velocity v\ of particle 1 at time t\ and the corresponding retarded 
and advanced velocities v 2b and v 2a of particle 2. 

Fig. 3 A symmetric orbit in arbitrary units. The symmetric orbit is completely defined 
by only two functions r(t) and v(t) , from which the quantities r a , r b , v±, v 2 b and v 2a are 
determined. Notice that — v 2a when r b = r is the same as v 1 when r a = r. 

Fig. 4 A symmetric orbit in arbitrary units. Illustrated are the segments of trajectory 
relevant to the evaluation of Y a and Y b and the quantities r + ,r~, v£ and v 2 ■ 

Fig. 5 Numerically calculated trajectories using the Pade approximation, in units where 
c = e = m = 1. For each value of the energy E, the regularized integrator starts from the 
initial condition r a — r b — 0, vi — 1, Y a — and Y b = . The four different orbits shown 
have energies E = —1.0, 0.1, 0.5 and 0.8 as indicated. 

Fig. 6. Numerically calculated 9{vi) (as defined by Eq. (54) ) plotted versus the velocity 
vi, in units where c = e = m = 1. This figure illustrates the jump in 9 near v — 1 at E — 0.8 
(Dash Dot), E = 0.9 (Dot) , E = 0.95 (Dash) and E = 0.99 (solid line). 
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Fig. 7. Numerically calculated 9{y\) (as defined by Eq. (54) ) plotted versus the velocity 
vi, in units where c — e — m — 1. In this figure the discontinuity region near v — 1 is 
blown-up. Energies are E = 0.8 (Dash Dot), E = 0.9 (Dot) , E = 0.95 (Dash) and E = 0.99 
(solid line). Notice that the jump in 6{v) at E = 0.99 is quite differently from the jump at 
the other energies. 

Fig. 8 Numerically calculated trajectories using the Spline interpolation, in units where 
c — e — m — 1. For each value of the energy E, the regularized integrator starts from the 
initial condition r a — — 0, V\ — 1, Y a — and Y b = \[2 . The four different orbits have 
energies E = 0.8, 0.9, 0.95 and 0.99 as indicated. 

Table 1: Numerically calculated orbits using the Pade approximation. Indicated are the 
energy E, the maximum light-cone distance r and the relative size of the first and last Pade 
coefficients, as well as the relative error of the minimization scheme, Ax/r Q . Notice that at 
the negative energy E = —1.0 the maximum radius r is less than the classical electronic 

2 

radius r* = = 1. 

Table 2: Numerically calculated orbits using the Spline interpolation. Indicated are 
the energy E, the maximum light-cone distance r , the saturation value r 9(v = 0) and 
the relative error of the minimization scheme, Ax/r^. Notice again that at E — —1.0 

2 

the maximum radius r is less than the classical electronic radius r* = = 1 and that 
r 9(v = 0) approximates the Coulombian limiting value of 3/2 as E tends to I. 
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